Household level determinants

Nguyen Bich Ngoc, Cedric Prevedello, Mario Cools, Jacques Teller 01 February 2022

## Warning in wtd.var(x, na.rm = T, us$weight): only one effective observation; variance
## estimate undefined

## Warning in wtd.var(x, na.rm = T, us$weight): only one effective observation; variance
## estimate undefined

## Warning in cor(itvar[, -(1:3)], use = "pairwise.complete.obs"): the standard
## deviation is zero

prop.table(table(us$rwusbn))

## 
##        no       yes 
## 0.5243039 0.4756961

prop.table(table(us$pvwlbn))

## 
##         no        yes 
## 0.95233601 0.04766399

## Warning: Removed 496 rows containing non-finite values (stat_boxplot).

## Warning: Removed 497 rows containing non-finite values (stat_boxplot).

## Warning: Removed 497 rows containing non-finite values (stat_boxplot).

## Warning: Removed 497 rows containing non-finite values (stat_smooth).

## Warning: Removed 497 rows containing missing values (geom_point).

## Warning: Removed 496 rows containing non-finite values (stat_smooth).

## Warning: Removed 496 rows containing missing values (geom_point).

## Warning: Removed 496 rows containing non-finite values (stat_boxplot).

## Warning: Removed 623 rows containing non-finite values (stat_boxplot).

mun_coord <- st_coordinates(st_centroid(mun_wal))

## Warning in st_centroid.sf(mun_wal): st_centroid assumes attributes are constant over
## geometries of x

knn <- knn2nb(knearneigh(mun_coord, k = 1), sym = T)
max_dist <- max(unlist(nbdists(knn, mun_coord)))
nb_dist <- dnearneigh(mun_coord, 0, (max_dist + 0.01), row.names = mun_wal$municd)
w_dist <- nb2listw(nb_dist, style = "W")

moran.test(x = mun_wal$csmptv, listw = w_dist, na.action = na.exclude, zero.policy = T)

## 
##  Moran I test under randomisation
## 
## data:  mun_wal$csmptv  
## weights: w_dist 
## omitted: 10, 79, 94, 96, 97, 121, 130, 131, 136, 144, 148, 168, 171, 173, 175, 176, 179, 180, 186, 195, 196, 197, 207, 215, 216, 219, 220, 223, 224, 226, 230, 231, 234, 238, 240, 241, 242, 243, 245, 247, 248, 249, 250, 253, 254, 255, 257, 259, 262   
## 
## Moran I statistic standard deviate = 3.089, p-value = 0.001004
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      0.0681561116     -0.0047169811      0.0005565335

cor.test(mun_wal$bltup2k5, mun_wal$rwt_prop, method = "spearman")

## Warning in cor.test.default(mun_wal$bltup2k5, mun_wal$rwt_prop, method = "spearman"):
## Cannot compute exact p-value with ties

## 
##  Spearman's rank correlation rho
## 
## data:  mun_wal$bltup2k5 and mun_wal$rwt_prop
## S = 2060316, p-value = 0.0001093
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.2614013

cor.test(mun_wal$bltup2k5, mun_wal$gard_prop, method = "spearman")

## Warning in cor.test.default(mun_wal$bltup2k5, mun_wal$gard_prop, method =
## "spearman"): Cannot compute exact p-value with ties

## 
##  Spearman's rank correlation rho
## 
## data:  mun_wal$bltup2k5 and mun_wal$gard_prop
## S = 2016113, p-value = 0.0005477
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.2343384

moran.test(x = mun_wal$rwt_prop, listw = w_dist, na.action = na.exclude, zero.policy = T)

## 
##  Moran I test under randomisation
## 
## data:  mun_wal$rwt_prop  
## weights: w_dist 
## omitted: 10, 79, 94, 96, 97, 121, 130, 131, 136, 144, 148, 168, 171, 173, 175, 176, 180, 186, 195, 196, 197, 207, 215, 216, 219, 220, 223, 224, 226, 230, 231, 234, 238, 240, 241, 242, 243, 245, 247, 248, 249, 250, 253, 254, 255, 257, 259, 262   
## 
## Moran I statistic standard deviate = 2.8048, p-value = 0.002517
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      0.0641176733     -0.0046948357      0.0006018964

moran.test(x = mun_wal$inceqa, listw = w_dist, na.action = na.exclude, zero.policy = T)

## 
##  Moran I test under randomisation
## 
## data:  mun_wal$inceqa  
## weights: w_dist 
## omitted: 10, 79, 94, 96, 97, 121, 130, 131, 136, 144, 148, 168, 171, 173, 175, 176, 180, 186, 195, 196, 197, 207, 215, 216, 219, 220, 223, 224, 226, 230, 231, 234, 238, 240, 241, 242, 243, 245, 247, 248, 249, 250, 253, 254, 255, 257, 259, 262   
## 
## Moran I statistic standard deviate = 1.8875, p-value = 0.02955
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      0.0414341086     -0.0046948357      0.0005972915

moran.test(x = mun_wal$hhs_tot, listw = w_dist, na.action = na.exclude, zero.policy = T)

## 
##  Moran I test under randomisation
## 
## data:  mun_wal$hhs_tot  
## weights: w_dist 
## omitted: 10, 79, 94, 96, 97, 121, 130, 131, 136, 144, 148, 168, 171, 173, 175, 176, 180, 186, 195, 196, 197, 207, 215, 216, 219, 220, 223, 224, 226, 230, 231, 234, 238, 240, 241, 242, 243, 245, 247, 248, 249, 250, 253, 254, 255, 257, 259, 262   
## 
## Moran I statistic standard deviate = 2.0369, p-value = 0.02083
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      0.0442581505     -0.0046948357      0.0005775723

cat("percentage of single member family")

## percentage of single member family

sum(us$hhs_tot %in% 1, na.rm = T)*100/sum(!is.na(us$hhs_tot))

## [1] 22.90192

cat("percentage of families with child(ren)")

## percentage of families with child(ren)

sum(us$hhs_0_17 > 0 & us$hhs_18_95 > 1, na.rm = T)*100/sum(!is.na(us$hhs_tot))

## [1] 22.14358

cat("percentage of couples without children")

## percentage of couples without children

sum(us$hhs_0_17 %in% 0 & us$hhs_18_95 %in% 2, na.rm = T)*100/sum(!is.na(us$hhs_tot))

## [1] 40.69767

cat("use distribution water for garden")

## use distribution water for garden

prop.table(table(us$garden))

## 
##       no      yes 
## 0.237848 0.762152

cat("presence of washing machine")

## presence of washing machine

prop.table(table(us$wasmch))

## 
##        no       yes 
## 0.0797546 0.9202454

cat("presence of dishwashing")

## presence of dishwashing

prop.table(table(us$dshwas))

## 
##       no      yes 
## 0.333176 0.666824

cat("(a) dual-flush or low-flush toilet(s)")

## (a) dual-flush or low-flush toilet(s)

prop.table(table(us$eftoil))

## 
##        no       yes 
## 0.2765455 0.7234545

cat("low flow shower head")

## low flow shower head

prop.table(table(us$efshhd))

## 
##        no       yes 
## 0.6083058 0.3916942

cat("confident")

## confident

prop.table(table(us$cfdiwq))

## 
##                        confident                 rather confident 
##                       0.50262030                       0.28156265 
## neither confident nor suspicious                rather suspicious 
##                       0.11434016                       0.04573606 
##                       suspicious                       no opinion 
##                       0.02382087                       0.03191996

cat("pay per volume")

## pay per volume

prop.table(table(us$ppusvl))

## 
##         no        yes 
## 0.01416431 0.98583569

# 4. baseline model ---------------

## 4.1. fitting -------------

### prepare data --------

full <-
  lm(
    csmptv ~ avrprc + hhs_18_95 + inceqa + parcarea + hhs_0_17 + nbbdrm + livapc + rwtuse + otsrid + rpjob + dshw3c + by_p1 + nbbtrm + dwowsh + pldisw + by_p2 + nbtoil + bthshw + buf300_max + rpgen + dfpay  + garden + prvncd + dwltyp + wasm3c + otsrod + rped4c + dwcsy9 + bldarea + efterp + eftech + nblvrm + ppusvl + dwelocp + btshrp + pvwlid + bdgmtr + hlpwtk + fseaid + nbktch + cfdiwq,
    data = fitdf,
    na.action = na.exclude,
    weights = weight
  )

anova(full)

## Analysis of Variance Table
## 
## Response: csmptv
##             Df Sum Sq Mean Sq  F value    Pr(>F)    
## avrprc       1  95594   95594  97.8221 < 2.2e-16 ***
## hhs_18_95    1 189215  189215 193.6247 < 2.2e-16 ***
## inceqa       1   2062    2062   2.1104 0.1466260    
## parcarea     1    804     804   0.8232 0.3644649    
## hhs_0_17     1 129405  129405 132.4208 < 2.2e-16 ***
## nbbdrm       3   4310    1437   1.4701 0.2211372    
## livapc       1      3       3   0.0034 0.9535901    
## rwtuse       3  94379   31460  32.1928 < 2.2e-16 ***
## otsrid       1   5668    5668   5.8005 0.0162053 *  
## rpjob       11  30145    2740   2.8043 0.0013167 ** 
## dshw3c       2   1630     815   0.8338 0.4346965    
## by_p1        1   1469    1469   1.5028 0.2205371    
## nbbtrm       3  11697    3899   3.9898 0.0077236 ** 
## dwowsh       3  25475    8492   8.6895 1.081e-05 ***
## pldisw       1   4299    4299   4.3991 0.0362136 *  
## by_p2        1    468     468   0.4786 0.4892275    
## nbtoil       3   3147    1049   1.0736 0.3593545    
## bthshw       3   9797    3266   3.3418 0.0187403 *  
## buf300_max   5  30480    6096   6.2381 1.049e-05 ***
## rpgen        1   1079    1079   1.1042 0.2936143    
## dfpay        1  12404   12404  12.6931 0.0003847 ***
## garden       1   4864    4864   4.9774 0.0259062 *  
## prvncd       1   1108    1108   1.1338 0.2872385    
## dwltyp       3   1813     604   0.6185 0.6030839    
## wasm3c       2    847     424   0.4336 0.6483093    
## otsrod       1  17062   17062  17.4597 3.197e-05 ***
## rped4c       3   2154     718   0.7348 0.5313609    
## dwcsy9       8  22073    2759   2.8235 0.0042274 ** 
## bldarea      1    262     262   0.2678 0.6049210    
## efterp       1    340     340   0.3476 0.5556144    
## eftech       1    622     622   0.6365 0.4251584    
## nblvrm       3  12494    4165   4.2617 0.0053102 ** 
## ppusvl       1   1337    1337   1.3684 0.2423802    
## dwelocp      1   5784    5784   5.9185 0.0151607 *  
## btshrp       3   6411    2137   2.1869 0.0880195 .  
## pvwlid       1    160     160   0.1641 0.6854587    
## bdgmtr       1    144     144   0.1475 0.7010445    
## hlpwtk       1     97      97   0.0993 0.7527233    
## fseaid       1   4729    4729   4.8392 0.0280521 *  
## nbktch       2   3811    1906   1.9500 0.1428199    
## cfdiwq       5  18833    3767   3.8543 0.0018306 ** 
## Residuals  979 956704     977                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1